Tau profiles

Tau is used here to calculate stage specificity

library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.3     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.4     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.0
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(myTAI)
specificity <- function(PES, measure = "tau", drop_Na = T){
        n_expcol <- ncol(PES) - 2 # e.g. stages or tissue
        ExpressionSet <- PES[,-1]
        if (measure == "tau") {
                norm_expression <- t(apply(
                        ExpressionSet[, -1], # geneID
                        1, 
                        function(row) row / max(row))) # normalise expression
                specificity_index <- data.frame(
                        # calculate Tau and return ExpressionSet with tau
                        tau = apply(
                                norm_expression,
                                1,
                                function(row) sum(1 - row) / (n_expcol-1)),
                        GeneID = ExpressionSet[, 1])
        }
        if (drop_Na) {
                specificity_index <- na.omit(specificity_index)
        }
        return(specificity_index)
}

Load PES

Non-transformed

Fd_PES <-
  readr::read_csv(file = "data/Fd_PES.csv") %>% select(-c(gamete,matSP))
## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES <-
  readr::read_csv(file = "data/Fs_PES.csv") %>% select(-c(gamete,matSP))
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

sqrt-tranformed

Fd_PES.sqrt <-
  readr::read_csv(file = "data/Fd_PES.sqrt.csv") %>% select(-c(gamete,matSP))
## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.sqrt <-
  readr::read_csv(file = "data/Fs_PES.sqrt.csv") %>% select(-c(gamete,matSP))
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

log2-tranformed

Fd_PES.log2 <-
  readr::read_csv(file = "data/Fd_PES.log2.csv") %>% select(-c(gamete,matSP))
## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.log2 <-
  readr::read_csv(file = "data/Fs_PES.log2.csv") %>% select(-c(gamete,matSP))
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

rank-tranformed

Fd_PES.rank <-
  readr::read_csv(file = "data/Fd_PES.rank.csv") %>% select(-c(gamete,matSP))
## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.rank <-
  readr::read_csv(file = "data/Fs_PES.rank.csv") %>% select(-c(gamete,matSP))
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

rlog-tranformed

Fd_PES.rlog <-
  readr::read_csv(file = "data/Fd_PES.rlog.csv") %>% select(-c(gamete,matSP))
## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.rlog <-
  readr::read_csv(file = "data/Fs_PES.rlog.csv") %>% select(-c(gamete,matSP))
## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

Get tau values

Fd_PES_tau <- specificity(Fd_PES)
Fs_PES_tau <- specificity(Fs_PES)
Fd_PES_tau.sqrt <- specificity(Fd_PES.sqrt)
Fs_PES_tau.sqrt <- specificity(Fs_PES.sqrt)
Fd_PES_tau.log2 <- specificity(Fd_PES.log2)
Fs_PES_tau.log2 <- specificity(Fs_PES.log2)
Fd_PES_tau.rank <- specificity(Fd_PES.rank)
Fs_PES_tau.rank <- specificity(Fs_PES.rank)

Plot tau distribution

hist(Fd_PES_tau$tau, breaks = 100)

hist(Fs_PES_tau$tau, breaks = 100)

hist(Fd_PES_tau.sqrt$tau, breaks = 100)

hist(Fs_PES_tau.sqrt$tau, breaks = 100)

hist(Fd_PES_tau.log2$tau, breaks = 100)

hist(Fs_PES_tau.log2$tau, breaks = 100)

hist(Fd_PES_tau.rank$tau, breaks = 100)

hist(Fs_PES_tau.rank$tau, breaks = 100)

Make tau into quantiles (deciles)

This is analogous to divergence map

tau_map <- function(tauExpressionSet, n_quantile = 10){
    GeneID <- tau_strata <- NULL
    data.table::setDT(tauExpressionSet)
    data.table::setkeyv(tauExpressionSet, c("tau", "GeneID"))
    tau_tbl_tauMap <- data.table::as.data.table(dplyr::select(dtplyr::lazy_dt(tauExpressionSet), 
        tau, GeneID))
    QuantileValues <- stats::quantile(tau_tbl_tauMap[, tau], 
        probs = seq(0, 1, (1/{
            {
                n_quantile
            }
        })), na.rm = TRUE)
    if (any(duplicated(QuantileValues))) {
        stop("ERROR: there is at least one duplicate threshold in the quantile analysis. Please select a lower value as n_quantile.")
    }
    else {
    }
    tau_tbl_tauMap$tau <- base::cut(tau_tbl_tauMap[, tau], 
        breaks = QuantileValues, include.lowest = TRUE, labels = FALSE)
    data.table::setnames(tau_tbl_tauMap, old = "tau", new = "tau_strata")
    return(tau_tbl_tauMap[, list(tau_strata, GeneID)])
    
}
Fd_PES_tau.quantile <- tau_map(Fd_PES_tau)
Fs_PES_tau.quantile <- tau_map(Fs_PES_tau)
Fd_PES_tau.sqrt.quantile <- tau_map(Fd_PES_tau.sqrt)
Fs_PES_tau.sqrt.quantile <- tau_map(Fs_PES_tau.sqrt)
Fd_PES_tau.log2.quantile <- tau_map(Fd_PES_tau.log2)
Fs_PES_tau.log2.quantile <- tau_map(Fs_PES_tau.log2)
Fd_PES_tau.rank.quantile <- tau_map(Fd_PES_tau.rank)
Fs_PES_tau.rank.quantile <- tau_map(Fs_PES_tau.rank)
dplyr::left_join(Fd_PES_tau, Fd_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fs_PES_tau, Fs_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fd_PES_tau.sqrt, Fd_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fs_PES_tau.sqrt, Fs_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fd_PES_tau.log2, Fd_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fs_PES_tau.log2, Fs_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

# box plots instead

dplyr::left_join(Fd_PES_tau, Fd_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(Fs_PES_tau, Fs_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(Fd_PES_tau.sqrt, Fd_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(Fs_PES_tau.sqrt, Fs_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(Fd_PES_tau.log2, Fd_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(Fs_PES_tau.log2, Fs_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(Fd_PES_tau.quantile, Fd_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau_strata)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fs_PES_tau.quantile, Fs_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau_strata)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fd_PES_tau.sqrt.quantile, Fd_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau_strata)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fs_PES_tau.sqrt.quantile, Fs_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau_strata)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fd_PES_tau.log2.quantile, Fd_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau_strata)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

dplyr::left_join(Fs_PES_tau.log2.quantile, Fs_PES) %>%
        ggplot2::ggplot(aes(x = PS, y = tau_strata)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_smooth(method = "lm")
## Joining with `by = join_by(GeneID)`
## `geom_smooth()` using formula = 'y ~ x'

corr_from_tauPES <- function(PES, PES_tau){
        input_data <- dplyr::left_join(
                PES_tau, 
                dplyr::select(PES, c(1, 2)), by = "GeneID")
        return(stats::cor.test(input_data$PS, input_data$tau_strata))
}
corr_from_tauPES(PES_tau = Fd_PES_tau.quantile, PES = Fd_PES)
## 
##  Pearson's product-moment correlation
## 
## data:  input_data$PS and input_data$tau_strata
## t = 5.0421, df = 7840, p-value = 4.708e-07
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.03476303 0.07888669
## sample estimates:
##        cor 
## 0.05685262
corr_from_tauPES(PES_tau = Fs_PES_tau.quantile, PES = Fs_PES)
## 
##  Pearson's product-moment correlation
## 
## data:  input_data$PS and input_data$tau_strata
## t = 3.9597, df = 8133, p-value = 7.569e-05
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.02215486 0.06553341
## sample estimates:
##        cor 
## 0.04386481
corr_from_tauPES(PES_tau = Fd_PES_tau.sqrt.quantile, PES = Fd_PES)
## 
##  Pearson's product-moment correlation
## 
## data:  input_data$PS and input_data$tau_strata
## t = 5.2219, df = 7840, p-value = 1.817e-07
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.03678741 0.08090072
## sample estimates:
##       cor 
## 0.0588728
corr_from_tauPES(PES_tau = Fs_PES_tau.sqrt.quantile, PES = Fs_PES)
## 
##  Pearson's product-moment correlation
## 
## data:  input_data$PS and input_data$tau_strata
## t = 4.3566, df = 8133, p-value = 1.337e-05
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.02654927 0.06991026
## sample estimates:
##       cor 
## 0.0482525
corr_from_tauPES(PES_tau = Fd_PES_tau.log2.quantile, PES = Fd_PES)
## 
##  Pearson's product-moment correlation
## 
## data:  input_data$PS and input_data$tau_strata
## t = 12.425, df = 7840, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.1171916 0.1606039
## sample estimates:
##       cor 
## 0.1389645
corr_from_tauPES(PES_tau = Fs_PES_tau.log2.quantile, PES = Fs_PES)
## 
##  Pearson's product-moment correlation
## 
## data:  input_data$PS and input_data$tau_strata
## t = 11.218, df = 8133, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.1019871 0.1447872
## sample estimates:
##       cor 
## 0.1234446

The correlation between PS and tau_strata is very weak. Therefore, we have a statistically independent measure.

# For F. distichus

test_data <- left_join(Fd_PES_tau.sqrt.quantile, Fd_PES_tau.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 117.33, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.9716791
test_data <- left_join(Fd_PES_tau.log2.quantile, Fd_PES_tau.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 93.753, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.7764369
test_data <- left_join(Fd_PES_tau.sqrt, Fd_PES_tau, by = "GeneID")
plot(test_data$tau.x, test_data$tau.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau.x, test_data$tau.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau.x and test_data$tau.y
## z = 127.73, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##      tau 
## 0.961914
test_data <- left_join(Fd_PES_tau.log2, Fd_PES_tau, by = "GeneID")
plot(test_data$tau.x, test_data$tau.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau.x, test_data$tau.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau.x and test_data$tau.y
## z = 97.842, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##      tau 
## 0.736815
# For F. serratus

test_data <- left_join(Fs_PES_tau.sqrt.quantile, Fs_PES_tau.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 119.15, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.9688347
test_data <- left_join(Fs_PES_tau.log2.quantile, Fs_PES_tau.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 94.345, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.7671409
test_data <- left_join(Fs_PES_tau.sqrt, Fs_PES_tau, by = "GeneID")
plot(test_data$tau.x, test_data$tau.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau.x, test_data$tau.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau.x and test_data$tau.y
## z = 129.65, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.9586059
test_data <- left_join(Fs_PES_tau.log2, Fs_PES_tau, by = "GeneID")
plot(test_data$tau.x, test_data$tau.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau.x, test_data$tau.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau.x and test_data$tau.y
## z = 98.2, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.7260861

deciling improves the correlation between

Fd_SpES <- myTAI::MatchMap(Fd_PES_tau.quantile, Fd_PES[-1])
Fs_SpES <- myTAI::MatchMap(Fs_PES_tau.quantile, Fs_PES[-1])
Fd_SpES.sqrt <- myTAI::MatchMap(Fd_PES_tau.sqrt.quantile, Fd_PES.sqrt[-1])
Fs_SpES.sqrt <- myTAI::MatchMap(Fs_PES_tau.sqrt.quantile, Fs_PES.sqrt[-1])
Fd_SpES.log2 <- myTAI::MatchMap(Fd_PES_tau.log2.quantile, Fd_PES.log2[-1])
Fs_SpES.log2 <- myTAI::MatchMap(Fs_PES_tau.log2.quantile, Fs_PES.log2[-1])

Plot TSpI

(Transcriptome Specificity Index). We should keep in mind that we are only retaining the embryo stages.

#TI stands for transcriptome index
TI.preplot <- function(ExpressionSet, permutations = 1000){
  std <- 
    myTAI::bootMatrix(
      ExpressionSet = ExpressionSet, 
      permutations  = permutations) %>% 
    apply(2, stats::sd)
  TI.out <- 
    myTAI::TAI(ExpressionSet) %>%
    tibble::as_tibble(rownames = "Stage", colnames = c("TAI", "PS")) %>% 
    dplyr::bind_cols(as_tibble(std)) %>%
    dplyr::rename(TAI = 2, sd = 3)
  return(TI.out)
}
get_TAI <- function(SpES, ordered_stages){
  TAI <- 
    TI.preplot(SpES)
  TAI_out <- 
    TAI
  TAI_out$Stage <- base::factor(TAI_out$Stage, ordered_stages)
  return(TAI_out)
}
get_TAI(Fd_SpES, ordered_stages = colnames(Fd_SpES)[-c(1,2)]) %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "TPM",
       y = "TSpI")
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

get_TAI(Fs_SpES, ordered_stages = colnames(Fs_SpES)[-c(1,2)]) %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "TPM",
       y = "TSpI")
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

# sqrt transformation

get_TAI(Fd_SpES.sqrt, ordered_stages = colnames(Fd_SpES)[-c(1,2)]) %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "sqrt(TPM)",
       y = "TSpI")
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

get_TAI(Fs_SpES.sqrt, ordered_stages = colnames(Fs_SpES)[-c(1,2)]) %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "sqrt(TPM)",
       y = "TSpI")
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

# log transformation

get_TAI(Fd_SpES.log2, ordered_stages = colnames(Fd_SpES)[-c(1,2)]) %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus distichus",
       subtitle = "log2(TPM+1)",
       y = "TSpI")
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

get_TAI(Fs_SpES.log2, ordered_stages = colnames(Fs_SpES)[-c(1,2)]) %>%
  ggplot2::ggplot(
    aes(
      y = TAI, 
      x = Stage,
      group = 1)) +
  ggplot2::geom_ribbon(
    aes(ymin = TAI - sd,
        ymax = TAI + sd
        ), 
    colour = "#00A087FF", 
    fill = "#00A087FF", 
    alpha = 0.1) +
  ggplot2::geom_line(
    size = 2,
    lineend = "round",
    colour = "#00A087FF") +
  theme_classic() +
  ggsci::scale_colour_npg() +
  labs(title = "Fucus serratus",
       subtitle = "log2(TPM+1)",
       y = "TSpI")
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## New names:
## • `value` -> `value...2`
## • `value` -> `value...3`

Plotting relative expression by Tau strata.

plot_relative_expression_tau <- function(SpES, average = "mean", plot = T){
  averaging_function <- base::match.fun(average)        

  PES_in <- 
    SpES %>%
    tidyr::pivot_longer(!c(tau_strata, GeneID), values_to = "Expression", names_to = "Stage") %>%
    dplyr::group_by(tau_strata,Stage) %>%
    dplyr::summarise(mean_expression = averaging_function(Expression),
                     n_genes = n())
  
  max_expression <- 
    PES_in %>%
    dplyr::group_by(tau_strata) %>%
    dplyr::summarise(max_expression = max(mean_expression),
                     min_expression = min(mean_expression))

  # normalise expression between 0 and 1
  PES_relative_expression <-
    left_join(PES_in, max_expression, by = "tau_strata") %>%
    dplyr::mutate(relative_expression = (mean_expression-min_expression)/(max_expression-min_expression))
  
  PES_relative_expression$Stage <- 
    factor(PES_relative_expression$Stage, levels = unique(colnames(SpES)[-c(1,2)]))
  
  # return if plotting is not desired.
  if(!plot){return(PES_relative_expression)}
  
  PES_relative_expression <-
    dplyr::mutate(PES_relative_expression, tau_strata_category = case_when(
      tau_strata < 6 ~ "low tau",
      TRUE~ "high tau"
    ))
  
  ggplot2::ggplot(PES_relative_expression,
                  aes(x = Stage, 
                      y = relative_expression, 
                      group = tau_strata, 
                      # linetype = tau_strata_representativeness,
                      # linewidth = log2(n_genes),
                      colour = factor(tau_strata))) +
  # ggplot2::geom_line(alpha = 0.5) +
  ggplot2::geom_line(linewidth = 2, linetype = "twodash") +
  ggplot2::facet_wrap(vars(tau_strata_category)) +
  ggplot2::labs(y = "relative expression", colour = "tau_strata")
}
plot_relative_expression_tau(Fs_SpES) + ggplot2::labs(title = "Fs", subtitle = "TPM")
## `summarise()` has grouped output by 'tau_strata'. You can override using the
## `.groups` argument.

plot_relative_expression_tau(Fd_SpES) + ggplot2::labs(title = "Fd", subtitle = "TPM")
## `summarise()` has grouped output by 'tau_strata'. You can override using the
## `.groups` argument.

plot_relative_expression_tau(Fs_SpES.sqrt) + ggplot2::labs(title = "Fs", subtitle = "sqrt(TPM)")
## `summarise()` has grouped output by 'tau_strata'. You can override using the
## `.groups` argument.

plot_relative_expression_tau(Fd_SpES.sqrt) + ggplot2::labs(title = "Fd", subtitle = "sqrt(TPM)")
## `summarise()` has grouped output by 'tau_strata'. You can override using the
## `.groups` argument.

plot_relative_expression_tau(Fs_SpES.log2) + ggplot2::labs(title = "Fs", subtitle = "log2(TPM)")
## `summarise()` has grouped output by 'tau_strata'. You can override using the
## `.groups` argument.

plot_relative_expression_tau(Fd_SpES.log2) + ggplot2::labs(title = "Fd", subtitle = "log2(TPM)")
## `summarise()` has grouped output by 'tau_strata'. You can override using the
## `.groups` argument.

The relative expression of genes with high tau (specific) decreases during the waist of the hourglass compared to genes with low tau (broad).

Performing statistical analysis for the significance of the hourglass-shape from tau.

# Fucus serratus
Fs_SpES_tS <- 
  tfStability(
    Fs_SpES,
    transforms = c("none", "sqrt", "log2", "rank"),
    permutations = 500
    )
Fs_SpES_tS_processed <-
  tibble::as_tibble(Fs_SpES_tS, rownames = "transformation") %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "FlatLineTest")

# Fucus distichus
Fd_SpES_tS <- 
  tfStability(
    Fd_SpES,
    transforms = c("none", "sqrt", "log2", "rank"),
    permutations = 500
    )
Fd_SpES_tS_processed <-
  tibble::as_tibble(Fd_SpES_tS, rownames = "transformation") %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "FlatLineTest")
Fs_SpES_tS_processed %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 1
    ) +
  ggplot2::labs(
    title = "Fucus serratus",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Fd_SpES_tS_processed %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 1
    ) +
  ggplot2::labs(
    title = "Fucus distichus",
    subtitle = "FlatLineTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

# Fucus serratus
Fs_SpES_tS <- 
  tfStability(
    Fs_SpES,
    TestStatistic = "ReductiveHourglassTest",
    transforms = c("none", "sqrt", "log2", "rank"),
    modules = list(early = 1:2, mid = 3:4, late = 5),
    permutations = 500
    )
Fs_SpES_tS_processed <-
  tibble::as_tibble(Fs_SpES_tS, rownames = "transformation") %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "ReductiveHourglassTest")

Fs_SpES_tS.alt <- 
  tfStability(
    Fs_SpES,
    TestStatistic = "ReductiveHourglassTest",
    transforms = c("none", "sqrt", "log2", "rank"),
    modules = list(early = 1:3, mid = 4, late = 5),
    permutations = 500
    )
Fs_SpES_tS_processed.alt <-
  tibble::as_tibble(Fs_SpES_tS.alt, rownames = "transformation") %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "ReductiveHourglassTest")

# Fucus distichus
Fd_SpES_tS <- 
  tfStability(
    Fd_SpES,
    TestStatistic = "ReductiveHourglassTest",
    transforms = c("none", "sqrt", "log2", "rank"),
    modules = list(early = 1:4, mid = 5, late = 6),
    permutations = 500
    )
Fd_SpES_tS_processed <-
  tibble::as_tibble(Fd_SpES_tS, rownames = "transformation") %>%
  dplyr::rename("pval" = value) %>%
  dplyr::mutate(test = "ReductiveHourglassTest")
Fs_SpES_tS_processed %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 1
    ) +
  ggplot2::labs(
    title = "Fucus serratus",
    subtitle = "ReductiveHourglassTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Fs_SpES_tS_processed.alt %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 1
    ) +
  ggplot2::labs(
    title = "Fucus serratus (only 3w as conserved)",
    subtitle = "ReductiveHourglassTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

Fd_SpES_tS_processed %>%
  ggplot2::ggplot(
    aes(
      y = -log10(pval), 
      x = transformation)) +
  ggplot2::geom_col(width = 0.05) +
  ggplot2::geom_point(size = 5) +
  ggplot2::geom_hline(
    yintercept = -log10(0.05), 
    colour = "#E64B35FF",
    linetype = 'dashed',
    size = 1
    ) +
  ggplot2::labs(
    title = "Fucus distichus",
    subtitle = "ReductiveHourglassTest",
    x = "Transformation"
  ) +
  ggplot2::coord_flip() +
  ggplot2::theme_classic()

I see that for Fucus serratus, we cannot use the same category for the most conserved stage (1w, 3w) for the tau data. Instead, there is a shift.

Tissues

Can we look at the out of testis hypothesis by looking at tau?

Load PES

Non-transformed

Fd_PES.adult <-
  readr::read_csv(file = "data/Fd_PES_matSP.csv")
## Rows: 7907 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.adultM <-
  readr::read_csv(file = "data/Fs_PES_M_matSP.csv")
## Rows: 8291 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.adultF <-
  readr::read_csv(file = "data/Fs_PES_F_matSP.csv")
## Rows: 8291 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

sqrt-tranformed

Fd_PES.adult.sqrt <-
  readr::read_csv(file = "data/Fd_PES_matSP.sqrt.csv")
## Rows: 7907 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.adultM.sqrt <-
  readr::read_csv(file = "data/Fs_PES_M_matSP.sqrt.csv")
## Rows: 8291 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.adultF.sqrt <-
  readr::read_csv(file = "data/Fs_PES_F_matSP.sqrt.csv")
## Rows: 8291 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

log2-tranformed

Fd_PES.adult.log2 <-
  readr::read_csv(file = "data/Fd_PES_matSP.log2.csv")
## Rows: 7907 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.adultM.log2 <-
  readr::read_csv(file = "data/Fs_PES_M_matSP.log2.csv")
## Rows: 8291 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_PES.adultF.log2 <-
  readr::read_csv(file = "data/Fs_PES_F_matSP.log2.csv")
## Rows: 8291 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (5): PS, holdfast, reptip, stipe, vegtip
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_PES_tau.adult.quantile <- 
        specificity(Fd_PES.adult) |> tau_map()
Fs_PES_tau.adultM.quantile <- 
        specificity(Fs_PES.adultM) |> tau_map()
Fs_PES_tau.adultF.quantile <- 
        specificity(Fs_PES.adultF) |> tau_map()
Fd_PES_tau.adult.sqrt.quantile <- 
        specificity(Fd_PES.adult.sqrt) |> tau_map()
Fs_PES_tau.adultM.sqrt.quantile <-
        specificity(Fs_PES.adultM.sqrt) |> tau_map()
Fs_PES_tau.adultF.sqrt.quantile <- 
        specificity(Fs_PES.adultF.sqrt) |> tau_map()
Fd_PES_tau.adult.log2.quantile <- 
        specificity(Fd_PES.adult.log2) |> tau_map()
Fs_PES_tau.adultM.log2.quantile <- 
        specificity(Fs_PES.adultM.log2) |> tau_map()
Fs_PES_tau.adultF.log2.quantile <- 
        specificity(Fs_PES.adultF.log2) |> tau_map()

I wondered if we have the same results as with stage-specificity tau * PS when we use tissue-specificity tau. It is basically the same.

dplyr::left_join(specificity(Fd_PES.adult), Fd_PES.adult) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fs_PES.adultM), Fs_PES.adultM) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fs_PES.adultF), Fs_PES.adultF) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fd_PES.adult.sqrt), Fd_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fs_PES.adultM.sqrt), Fs_PES.adultM) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fs_PES.adultF.sqrt), Fs_PES.adultF) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fd_PES.adult.log2), Fd_PES) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fs_PES.adultM.log2), Fs_PES.adultM) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

dplyr::left_join(specificity(Fs_PES.adultF.log2), Fs_PES.adultF) %>%
        ggplot2::ggplot(aes(x = as.factor(PS), y = tau)) +
        ggplot2::geom_jitter(alpha = 0.1) +
        ggplot2::geom_boxplot(width = 0.2)
## Joining with `by = join_by(GeneID)`

Tau from adult tissue vs. from development - not strongly correlated.

# For F. distichus

test_data <- left_join(Fd_PES_tau.sqrt.quantile, Fd_PES_tau.adult.sqrt.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 44.507, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.3686045
test_data <- left_join(Fd_PES_tau.log2.quantile, Fd_PES_tau.adult.log2.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 56.334, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.4665532
# For F. serratus

test_data <- left_join(Fs_PES_tau.sqrt.quantile, Fs_PES_tau.adultM.sqrt.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 42.142, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.3433586
test_data <- left_join(Fs_PES_tau.log2.quantile, Fs_PES_tau.adultM.log2.quantile, by = "GeneID")
plot(test_data$tau_strata.x, test_data$tau_strata.y, col=rgb(0, 0, 0, 0.01))

stats::cor.test(test_data$tau_strata.x, test_data$tau_strata.y, method = "kendall")
## 
##  Kendall's rank correlation tau
## 
## data:  test_data$tau_strata.x and test_data$tau_strata.y
## z = 56.381, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.4593697

Plotting embryo stages with adult tissue specificity (to avoid statistical non-independence). The resulting effect is less pronounced but the pattern resembles an hourglass or late conservation. Of course, we are lacking resolution of conditions and it is unclear what one means when

myTAI::MatchMap(Fd_PES_tau.adult.quantile, Fd_PES[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.102  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultM.quantile, Fs_PES[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.101  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultF.quantile, Fs_PES[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.099  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fd_PES_tau.adult.sqrt.quantile, Fd_PES.sqrt[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.152  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultM.sqrt.quantile, Fs_PES.sqrt[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.103  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultF.sqrt.quantile, Fs_PES.sqrt[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.106  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fd_PES_tau.adult.log2.quantile, Fd_PES.log2[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.115  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultM.log2.quantile, Fs_PES.log2[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.119  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultF.log2.quantile, Fs_PES.log2[-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.189  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

Not shown but it looks slightly better

# myTAI::MatchMap(specificity(Fd_PES.adult), Fd_PES[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fs_PES.adultM), Fs_PES[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fs_PES.adultF), Fs_PES[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fd_PES.adult.sqrt), Fd_PES.sqrt[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fs_PES.adultM.sqrt) , Fs_PES.sqrt[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fs_PES.adultF.sqrt) , Fs_PES.sqrt[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fd_PES.adult.log2) , Fd_PES.log2[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fs_PES.adultM.log2) , Fs_PES.log2[-1]) |> myTAI::TAI() |> plot()
# myTAI::MatchMap(specificity(Fs_PES.adultM.log2) , Fs_PES.log2[-1]) |> myTAI::TAI() |> plot()
myTAI::MatchMap(Fd_PES_tau.adult.quantile, Fd_PES.adult[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.104  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  not significant (= no evolutionary signature in the transcriptome).
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultM.quantile, Fs_PES.adultM[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.152  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultF.quantile, Fs_PES.adultF[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.123  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fd_PES_tau.adult.sqrt.quantile, Fd_PES.adult.sqrt[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.139  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultM.sqrt.quantile, Fs_PES.adultM.sqrt[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.209  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultF.sqrt.quantile, Fs_PES.adultF.sqrt[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.114  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fd_PES_tau.adult.log2.quantile, Fd_PES.adult.log2[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.111  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultM.log2.quantile, Fs_PES.adultM.log2[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.119  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

myTAI::MatchMap(Fs_PES_tau.adultF.log2.quantile, Fs_PES.adultF.log2[,-1]) |> myTAI::PlotSignature(ylab = "TSpI")
## Plot signature: ' TAI ' and test statistic: ' FlatLineTest ' running  1000  permutations.
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## [ Number of Eigen threads that are employed on your machine: 8 ]
## 
## [ Computing age assignment permutations for test statistic ... ]
## 
[=========================================] 100%   
## [ Computing variances of permuted transcriptome signatures ... ]
## 
## 
## Total runtime of your permutation test: 0.119  seconds.
## 
## -> We recommended using at least 20000 permutations to achieve a sufficient permutation test.
## 
## Significance status of signature:  significant.
## 
## -> Now run 'FlatLineTest(..., permutations  = 1000, plotHistogram = TRUE)' to analyse the permutation test performance.

Get session info.

devtools::session_info()
## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value
##  version  R version 4.2.2 (2022-10-31)
##  os       macOS Big Sur ... 10.16
##  system   x86_64, darwin17.0
##  ui       X11
##  language (EN)
##  collate  en_US.UTF-8
##  ctype    en_US.UTF-8
##  tz       Europe/Berlin
##  date     2024-02-13
##  pandoc   3.1.6.2 @ /usr/local/bin/ (via rmarkdown)
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  package      * version    date (UTC) lib source
##  bit            4.0.5      2022-11-15 [1] CRAN (R 4.2.0)
##  bit64          4.0.5      2020-08-30 [1] CRAN (R 4.2.0)
##  bslib          0.5.1      2023-08-11 [1] CRAN (R 4.2.0)
##  cachem         1.0.8      2023-05-01 [1] CRAN (R 4.2.0)
##  callr          3.7.3      2022-11-02 [1] CRAN (R 4.2.0)
##  cli            3.6.1      2023-03-23 [1] CRAN (R 4.2.0)
##  codetools      0.2-19     2023-02-01 [1] CRAN (R 4.2.0)
##  colorspace     2.1-0      2023-01-23 [1] CRAN (R 4.2.0)
##  crayon         1.5.2      2022-09-29 [1] CRAN (R 4.2.0)
##  data.table     1.14.8     2023-02-17 [1] CRAN (R 4.2.0)
##  devtools       2.4.5      2022-10-11 [1] CRAN (R 4.2.0)
##  digest         0.6.33     2023-07-07 [1] CRAN (R 4.2.0)
##  dplyr        * 1.1.3      2023-09-03 [1] CRAN (R 4.2.0)
##  dtplyr         1.3.1      2023-03-22 [1] CRAN (R 4.2.0)
##  ellipsis       0.3.2      2021-04-29 [1] CRAN (R 4.2.0)
##  evaluate       0.22       2023-09-29 [1] CRAN (R 4.2.2)
##  fansi          1.0.5      2023-10-08 [1] CRAN (R 4.2.2)
##  farver         2.1.1      2022-07-06 [1] CRAN (R 4.2.0)
##  fastmap        1.1.1      2023-02-24 [1] CRAN (R 4.2.0)
##  fitdistrplus   1.1-11     2023-04-25 [1] CRAN (R 4.2.0)
##  forcats      * 1.0.0      2023-01-29 [1] CRAN (R 4.2.0)
##  foreach        1.5.2      2022-02-02 [1] CRAN (R 4.2.0)
##  fs             1.6.3      2023-07-20 [1] CRAN (R 4.2.0)
##  generics       0.1.3      2022-07-05 [1] CRAN (R 4.2.0)
##  ggplot2      * 3.4.4      2023-10-12 [1] CRAN (R 4.2.2)
##  ggsci          3.0.0      2023-03-08 [1] CRAN (R 4.2.0)
##  glue           1.6.2      2022-02-24 [1] CRAN (R 4.2.0)
##  gtable         0.3.4      2023-08-21 [1] CRAN (R 4.2.0)
##  hms            1.1.3      2023-03-21 [1] CRAN (R 4.2.0)
##  htmltools      0.5.6.1    2023-10-06 [1] CRAN (R 4.2.2)
##  htmlwidgets    1.6.2      2023-03-17 [1] CRAN (R 4.2.0)
##  httpuv         1.6.11     2023-05-11 [1] CRAN (R 4.2.2)
##  iterators      1.0.14     2022-02-05 [1] CRAN (R 4.2.0)
##  jquerylib      0.1.4      2021-04-26 [1] CRAN (R 4.2.0)
##  jsonlite       1.8.7      2023-06-29 [1] CRAN (R 4.2.0)
##  knitr          1.44       2023-09-11 [1] CRAN (R 4.2.0)
##  labeling       0.4.3      2023-08-29 [1] CRAN (R 4.2.0)
##  later          1.3.1      2023-05-02 [1] CRAN (R 4.2.2)
##  lattice        0.21-9     2023-10-01 [1] CRAN (R 4.2.2)
##  lifecycle      1.0.3      2022-10-07 [1] CRAN (R 4.2.0)
##  lubridate    * 1.9.3      2023-09-27 [1] CRAN (R 4.2.0)
##  magrittr       2.0.3      2022-03-30 [1] CRAN (R 4.2.0)
##  MASS           7.3-60     2023-05-04 [1] CRAN (R 4.2.2)
##  Matrix         1.5-4.1    2023-05-18 [1] CRAN (R 4.2.0)
##  memoise        2.0.1      2021-11-26 [1] CRAN (R 4.2.0)
##  mgcv           1.9-0      2023-07-11 [1] CRAN (R 4.2.0)
##  mime           0.12       2021-09-28 [1] CRAN (R 4.2.0)
##  miniUI         0.1.1.1    2018-05-18 [1] CRAN (R 4.2.0)
##  munsell        0.5.0      2018-06-12 [1] CRAN (R 4.2.0)
##  myTAI        * 1.0.1.9000 2023-12-07 [1] Github (drostlab/myTAI@27a2639)
##  nlme           3.1-163    2023-08-09 [1] CRAN (R 4.2.0)
##  pillar         1.9.0      2023-03-22 [1] CRAN (R 4.2.0)
##  pkgbuild       1.4.2      2023-06-26 [1] CRAN (R 4.2.0)
##  pkgconfig      2.0.3      2019-09-22 [1] CRAN (R 4.2.0)
##  pkgload        1.3.3      2023-09-22 [1] CRAN (R 4.2.0)
##  prettyunits    1.2.0      2023-09-24 [1] CRAN (R 4.2.0)
##  processx       3.8.2      2023-06-30 [1] CRAN (R 4.2.0)
##  profvis        0.3.8      2023-05-02 [1] CRAN (R 4.2.0)
##  promises       1.2.1      2023-08-10 [1] CRAN (R 4.2.2)
##  ps             1.7.5      2023-04-18 [1] CRAN (R 4.2.0)
##  purrr        * 1.0.2      2023-08-10 [1] CRAN (R 4.2.2)
##  R6             2.5.1      2021-08-19 [1] CRAN (R 4.2.0)
##  Rcpp           1.0.11     2023-07-06 [1] CRAN (R 4.2.0)
##  readr        * 2.1.4      2023-02-10 [1] CRAN (R 4.2.0)
##  remotes        2.4.2.1    2023-07-18 [1] CRAN (R 4.2.2)
##  rlang          1.1.1      2023-04-28 [1] CRAN (R 4.2.0)
##  rmarkdown      2.25       2023-09-18 [1] CRAN (R 4.2.2)
##  rstudioapi     0.15.0     2023-07-07 [1] CRAN (R 4.2.0)
##  sass           0.4.7      2023-07-15 [1] CRAN (R 4.2.0)
##  scales         1.2.1      2022-08-20 [1] CRAN (R 4.2.0)
##  sessioninfo    1.2.2      2021-12-06 [1] CRAN (R 4.2.0)
##  shiny          1.7.5.1    2023-10-14 [1] CRAN (R 4.2.2)
##  stringi        1.7.12     2023-01-11 [1] CRAN (R 4.2.0)
##  stringr      * 1.5.0      2022-12-02 [1] CRAN (R 4.2.0)
##  survival       3.5-7      2023-08-14 [1] CRAN (R 4.2.0)
##  tibble       * 3.2.1      2023-03-20 [1] CRAN (R 4.2.0)
##  tidyr        * 1.3.0      2023-01-24 [1] CRAN (R 4.2.0)
##  tidyselect     1.2.0      2022-10-10 [1] CRAN (R 4.2.0)
##  tidyverse    * 2.0.0      2023-02-22 [1] CRAN (R 4.2.0)
##  timechange     0.2.0      2023-01-11 [1] CRAN (R 4.2.0)
##  tzdb           0.4.0      2023-05-12 [1] CRAN (R 4.2.2)
##  urlchecker     1.0.1      2021-11-30 [1] CRAN (R 4.2.0)
##  usethis        2.2.2      2023-07-06 [1] CRAN (R 4.2.0)
##  utf8           1.2.3      2023-01-31 [1] CRAN (R 4.2.0)
##  vctrs          0.6.4      2023-10-12 [1] CRAN (R 4.2.2)
##  vroom          1.6.4      2023-10-02 [1] CRAN (R 4.2.2)
##  withr          2.5.1      2023-09-26 [1] CRAN (R 4.2.0)
##  xfun           0.40       2023-08-09 [1] CRAN (R 4.2.2)
##  xtable         1.8-4      2019-04-21 [1] CRAN (R 4.2.0)
##  yaml           2.3.7      2023-01-23 [1] CRAN (R 4.2.0)
## 
##  [1] /Library/Frameworks/R.framework/Versions/4.2/Resources/library
## 
## ──────────────────────────────────────────────────────────────────────────────